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1.  Introduction 


Instrumented  indentation  techniques  enable  characterization  of  the  mechanical  response  of  small 
material  samples  over  wide  ranges  of  spatial  resolution  and  deformation  magnitude  via  choice  of 
size  and  shape  of  indenters  and  depth  of  indentation.  Indentation  methods  have  been  used  to 
study,  for  example,  dislocation  accumulation  and  strain  gradient  effects  in  ductile  metals  (Nix  & 
Gao,  1998),  lattice  friction  resistance  to  dislocation  mobility  in  body-centered  cubic  crystals  (Qiu 
et  ah,  2001),  fracture  processes  in  ceramics  (Clarke  &  Tandon,  1995),  and  the  mechanical 
behavior  of  thin  films  (Saha  &  Nix,  2002;  Gerberich  et  ah,  2003).  The  underlying  theory  used  to 
derive  mechanical  properties  from  indentation  data  (Oliver  &  Pharr,  1992;  Hay  et  ah,  1999)  has 
historically  been  based  upon  isothermal,  rate-independent  elastic  or  elastic-plastic  constitutive 
assumptions — the  latter  extending  to  the  inelastic  response  of  the  original  solution  of  Hertz  (1896) 
of  elastic  contact  between  parabolic  solids  of  revolution.  In  recent  years,  VanLandingham  and 
coworkers  (VanLandingham  et  ah,  2001;  VanLandingham,  2003)  have  developed  a  dynamic 
indentation  method  exploiting  frequency  response  characteristics  of  the  substrate  material  in  order 
to  probe  time-dependent  constitutive  behavior  at  large  time  scales  (e.g.,  viscoelasticity  and  creep). 
This  method  involves  a  superimposition  of  high  frequency,  small  amplitude  loadings  with  a  single 
quasi-static  cycle  of  larger  magnitude.  Experimental  techniques  exploiting  moire  interferometry 
(Lu  et  al.,  2003)  have  likewise  been  developed  for  elucidating  viscous  effects  at  small  time  scales 
(e.g.,  high  rate  dynamic  plasticity). 

Numerical  simulations  based  on  finite  element  analysis  provide  opportune  support  for  the 
process  of  development  of  physically  realistic  constitutive  models,  and  the  acquisition  of 
necessary  material  parameters  for  these  models  via  matching  of  simulated  indentation  responses 
with  force  displacement  data  acquired  from  physical  experiments  (Bhattacharya  &  Nix,  1988). 
Typically,  numerical  indentation  studies  are  conducted  to  acquire  elastic  moduli,  hardness 
values,  and  less  commonly  in  plastic  materials,  the  static  yield  point  and  strain  hardening 
exponent  (Dao  et  al.,  2001;  Chollacoop  et  al.,  2003).  Hill  and  coworkers  (Hill  et  al.,  1989) 
performed  a  detailed  theoretical  and  numerical  study  of  the  Brinell  hardness  test  (i.e.,  indentation 
by  a  spherical  ball)  applied  to  elastic-plastic  media.  Vaidyanathan  et  al.  (2001)  used  a  combined 
experimental-numerical  approach  to  characterize  pressure-sensitive  yield  in  a  metallic  glass. 
Methods  have  also  been  developed  for  reproducing  entire  uniaxial  stress-strain  curves  purely 
from  spherical  indentation  load  displacement  data  (Herbert  et  al.,  2001),  although  these  are 
presently  regarded  by  the  scientific  community  as  rough  estimates  at  best. 

In  the  present  work,  we  study  the  response  of  elastic-plastic  materials  subjected  to  spherical 
indentation.  Section  2  describes  the  problem  setup,  specifically  an  axisymmetric  finite  element 
model  of  a  rigid  indenter  contacting  a  thick,  homogeneous  slab  of  substrate  material.  A  finite 
strain,  rate-independent,  isotropic  elastic-plastic  constitutive  formulation  is  invoked  for  the 
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behavior  of  the  substrate.  In  section  3,  the  relative  effects  of  selected  simulation  parameters  are 
investigated,  including  mesh  geometry,  elastic  stiffness  constants,  and  plastic  properties.  A 
titanium  alloy  (Ti-6A1-4V)  was  chosen  as  the  default  substrate  material  in  this  parametric 
investigation.  In  section  4,  we  compare  results  of  numerical  simulations  with  experimental  load 
displacement  data  for  alloyed  and  polymeric  blends  comprised  of  tungsten  (W),  tin  (Sn),  and/or 
nylon.  Conclusions  and  recommendations  follow  in  section  5. 


2.  Model 


2.1  Constitutive  Formulation 

The  finite  element  model  consists  of  a  rigid  spherical  indenter  and  elastoplastic  substrate.  The 
mechanical  response  of  the  substrate  is  dictated  by  classical,  rate-independent,  isotropic 
elastoplasticity  theory.  Only  minimal  aspects  of  this  theory  are  given  here;  for  additional  details, 
the  reader  may  consult  the  literature  (Hill,  1950;  Marin  &  McDowell,  1996).  Hereafter,  vectors 
and  tensors  of  rank  2  and  above  are  written  in  bold  face  type,  and  scalars  are  written  in  italics. 

The  kinematic  description  is  based  on  an  additive  decomposition  of  the  velocity  gradient  L,  an 
assumption  consistent  with  the  multiplicative  decomposition  of  Lee  (1969)  only  when  the 
condition  of  small  elastic  stretch  is  applicable: 

dv 

L  =  —  =  D  +  W  =  De+Dp+W,  (1) 

dx 

in  which  x  are  spatial  coordinates,  v  =  x  with  the  superposed  dot  a  material  time  derivative,  the 
spin  2W  =  L  -  Lr  with  the  superposed  T  denoting  the  transposition,  and  De  and  Dp  the  elastic 
and  plastic  deformation  rates,  respectively.  Following  the  introduction  of  suitable  co-rotational 
measures,  the  strain  rate  decomposition  may  be  written  in  time-integrated  form  via  the 
transformation 


D  =  De  +  Dp  — >  £  =  £e  +  £p  ,  (2) 

where  the  elastic  (ee)  and  plastic  (8P)  strain  tensors  sum  to  the  total  integrated  strain  measure  8, 
the  latter  equivalent  to  the  symmetric  part  of  the  spatial  displacement  gradient  only  in  the  limit  of 
small  strains  and  rotations.  The  isotropic  elastic  response  is  described  by 

o/  =  2G£e/,  P  =  ~K£v,  (3) 

where  cf  =  o  +  pi  is  the  Cauchy  stress  deviator  ,p  =  -(l/3)/r(o)is  the  hydrostatic  pressure, 

£e/  =  £e  -  (1  /  3)f  l ,  and  £v  =  tr(z) ,  assuming  isochoric  plastic  deformation.  The  identity  tensor 


2 


is  written  here  as  1.  Shear  modulus  G,  bulk  modulus  K,  Young’s  modulus  E,  and  Poisson’s  ratio 
V  are  interrelated  via 

2G(\  +  v)  =  3K(l-2v)  =  E .  (4) 


The  flow  rule  for  plastic  strain  increment  d£p  over  time  sub-step  At  is  written 

d£v  =  epnpAt , 


(5) 


in  which  ep  is  the  scalar  effective  plastic  strain  rate  and  the  direction  of  flow  is  coaxial  with  the 
deviatoric  stress,  i.e., 


n 


p 


(6) 


with  the  effective  stress  2a2  =  3a' :  a' .  The  plastic  strain  rate  is  determined  by  consistency  with 
the  yield  condition 


a  =  ay(ep), 


(V) 


where  the  yield  stress  <jv  is  chosen  here  as  a  simple  monotonic  function  of  the  cumulative 
effective  plastic  strain  ep  =  J  epdt . 


2.2  Problem  Setup 

Finite  element  analyses  were  conducted  with  the  ABAQUS '/Implicit  (2002)  commercial  package, 
with  the  preceding  constitutive  formulation  available  in  the  software’s  standard  material  library. 
All  analyses  were  isothermal  and  quasi-static  (i.e.,  inertial  forces  were  neglected).  As  shown  in 
figure  1,  the  boundary  value  problem  exhibits  an  axisymmetric  geometry  in  the  r  -  z  plane.  The 
indenter  radius  is  labeled  R,  and  W  and  H  denote  the  substrate  dimensions.  The  initial  configura¬ 
tion  is  shown  in  figure  1(a).  During  the  simulation  history,  the  indenter  is  depressed  in  the  -z- 
direction  a  distance  d  into  the  substrate  (figure  1(b)),  and  then  returned  along  the  reverse  path  to 
its  initial  position  (figure  1(c)).  The  force  acting  on  the  indenter,  directed  in  the  -z-direction,  is 
labeled  P. 


R 

1 — H 

w 

P 

1 

K 

¥ 

H 

i 

d 

W 

± 

d 

r 

(a)  (b)  (c) 


Figure  1.  Initial  (a),  intermediate  (b),  and  final  (c)  indenter  positions. 


1 ABAQUS  is  not  an  acronym. 
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The  indenter  was  modeled  as  an  analytical  rigid  surface,  while  the  deformable  substrate  was 
discretized  into  four-node  quadrilateral  (axisymmetric)  elements  featuring  reduced  integration 
and  “hourglass”  control.  Boundary  conditions  on  the  substrate  were  specified  by 

ur  =uz  =  0  (along  z  =  0) 

u  =  0  (along  r  =  0) 

r  =  <T rznr  =  0  (along  r  =  W)  (8) 

tz  =  azzn .  =  Ol 

>  (along  contact-free  portion  of  z  =  H ) 
tr  =  <Jrznz  =  Oj 

in  which  u  is  the  displacement  field,  t  is  the  traction  vector  per  current  configuration  area,  and  n 
is  the  unit  normal  vector  to  the  surface  of  the  body.  Preliminary  simulations  were  conducted  in 
order  to  investigate  effects  of  indenter  and  substrate  dimensions,  finite  element  mesh  densities, 
contact  friction,  and  elastoplastic  properties,  as  discussed  in  the  next  section.  Various  boundary 
conditions  were  investigated  for  the  far-field  edge  defined  by  r=W,  also  discussed  in  section  3. 

A  typical  mesh  is  shown  in  figure  2,  here  consisting  of  4900  elements.  Notice  that  the  mesh 
density  is  deliberately  highest  in  the  vicinity  of  the  contact  region  beneath  the  tip  of  indenter. 


Figure  2.  Finite  element  mesh,  4900  elements. 
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3.  Preliminary  Simulations 


The  default  material  of  study  in  the  preliminary  simulations  was  Ti-6A1-4V,  a  ductile  alloy 
featuring  a  polycrystalline  hexagonal  closely  packed  (HCP)  microstructure.  The  elastic 
constants  E  and  V  and  effective  deviatoric  yield  stress  0  =  0 r  (ep)  for  this  material  are  given  in 

figure  3,  obtained,  respectively,  from  wave  speed  and  mass  density  measurements  and  a  uniaxial 
quasi-static  tension  test,  with  £  the  logarithmic  strain  (see  Burkins  et  ah,  2001).  Note  that  o  is 
equivalent  to  the  measured  true  stress  in  a  uniaxial  stress  state.  In  the  simulations,  the  elastic 
constants  and  post-yielding  stress-strain  behavior  were  specified  as  separate  input  data  (see 
appendix  A). 


Figure  3.  Elastic  and  plastic  mechanical  properties  for 
Ti-6A1-4V. 


Figure  4  shows  the  effective  deviatoric  stress  in  the  substrate,  in  the  loaded  and  unloaded 
configurations  corresponding  to  figure  1(b)  and  1(c),  respectively.  Notice  that  while  the  stress 
field  is  supported  by  several  hundred  elements,  the  response  does  not  spread  to  the  lower  and 
right-most  edges  of  the  mesh,  along  z  =  0  and  r=W,  respectively.  Contours  in  figure  4(a) 
correspond  to  the  peak  load  applied  to  the  indenter,  while  those  in  figure  4(b)  denote  residual 
stresses  remaining  in  the  substrate  upon  full  contact  release. 

Figure  5  depicts  the  total  force  acting  on  the  rigid  indenter  in  the  -z  -direction,  P,  versus  the  total 
displacement  of  the  indenter  in  the  -z-direction,  8,  for  various  mesh  densities  and  far- (I eld 
boundary  conditions.  In  the  simulations,  c>  — >  d  during  the  loading  phase,  and  then  c>  — >  0  during 
unloading.  Note  that  such  P  versus  8 data  is  our  key  result  acquired  from  the  simulations,  as  it  is 
used  extensively  later  to  compare  experimental  and  simulated  findings.  In  figure  5,  “fixed  edge” 
denotes  the  displacement  condition  u(r  =  W)  =  0,  while  “free  edge”  denotes  the  traction  condition 
t'(r  =  W)  =  0.  From  figure  5,  we  conclude  that  the  choice  of  far-field  boundary  condition  has  a 
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negligible  effect  on  the  measured  P-8  response;  thus,  the  condition  u(r  =  W)  =  0  was  used 
exclusively  in  all  other  calculations  reported  here.  Analogously,  specification  of  u  =  0  along  the 
base  of  the  substrate  (along  z  =  0)  was  deemed  appropriate  (as  opposed  to  a  null  shear  stress 
boundary  condition),  as  the  stress  field  did  not  interact  with  this  boundary  in  the  present 
simulations.  Note  that  approximately  20  elements  at  the  surface  are  in  direct  contact  with  the 
indenter  in  figure  4(a).  Furthermore,  while  2500  elements  appear  to  offer  sufficient  resolution, 
the  4900-element  mesh  was  used  in  all  subsequent  calculations,  taking  advantage  of  short 
execution  times  for  the  analyses  (<900  seconds  on  a  single  R14000  Silicon  Graphics, 
Incorporated,  processor  for  the  4900-element  mesh  upon  indentation  to  d  =  5  //m). 
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Figure  4.  Stress  <7  [Pa]  in  (a)  deformed  intermediate  and  (b)  unloaded  final  configurations,  R  =  25  //m,  d  =  2  //m 
(frictionless  contact,  with  u\r  =  W]  =  0,  4900-element  mesh). 


Figure  5.  indentation  force  P  versus  displacement  S.  effects  of  mesh 
densities  and  far-field  boundary  conditions. 
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Figure  6  illustrates  effects  of  indenter  radius  R  and  coulomb  friction  coefficient  p  =  zip,  in  which 
T  is  the  resultant  shear  stress  opposing  relative  tangential  motion  of  two  surfaces  sharing  contact 
pressure  p.  The  geometry  of  the  indenter  clearly  imposes  a  dramatic  effect  on  the  P-8  response, 
which  is  an  important  consideration  when  one  is  matching  simulated  and  experimental  data, 
given  the  small  size  of  the  experimental  apparatus  whose  tip  shape  and  size  may  not  conform 
precisely  to  advertised  machining  specifications.  Friction  is  apparently  less  important,  at  least  in 
the  context  of  our  simulations,  which  assume  an  infinitely  rigid  indenter.  For  this  reason,  and 
since  frictionless  contact  is  the  norm  in  other  numerical  studies  reported  in  the  literature 
(Bhattacharya  &  Nix,  1988;  Dao  et  ah,  2001;  Vaidyanathan  et  ah,  2001;  Chollacoop  et  ah,  2003), 
p  =  0  is  invoked  in  all  remaining  simulations  of  the  present  work. 


Figure  6.  Indentation  force  P  versus  displacement  S. 

effects  of  indenter  radius  and  friction  coefficient, 
d  =  2  pm. 

Figures  7  and  8  delineate  the  effects  of  elastic  constants  E  and  v,  respectively,  on  the  computed 
P-8  response.  For  simulations  corresponding  to  figures  7(a)  and  8(a),  the  final  indentation  depth 
d  =  5  pm  instills  significant  yielding  and  finite  plastic  deformation  in  the  indentation  region, 
while  for  simulations  corresponding  to  figures  7(b)  and  8(b),  the  relatively  shallow  final  depth  of 
d  =  0.05  pm  induces  little  or  no  plastic  flow,  with  the  presence  and  degree  of  yielding  depending 
on  the  choice  of  elastic  constants.  In  the  simulations  associated  with  figures  7  and  8,  all  material 
parameters  (e.g.,  flow  stress  versus  plastic  strain)  corresponded  to  those  of  the  actual  Ti-6A1-4V 
material  of  figure  3,  except  for  the  elastic  constants  that  were  varied  parametrically  as  indicated 
in  the  figure  legends  (in  the  legends,  oy  =  <TV  \eP  =  0]  =  1000  MPa).  From  figure  7,  we  notice  the 
strong  effect  of  E  on  the  unloading  portion  of  each  P- <7 curve.  From  figure  8,  the  indentation 
load  increases  slightly  in  conjunction  with  decreasing  compressibility  of  the  material  (i.e.,  with 
increasing  v).  When  we  compare  figures  7  and  8,  Poisson’s  ratio  appears  to  exert  a  lesser  effect 
on  the  response  than  Young’s  modulus,  in  agreement  with  analytical  theory  (Oliver  &  Pharr, 
1992).  Furthermore,  the  effects  of  elastic  constants  decrease  with  increasing  indentation  depth, 
as  inelastic  mechanisms  increasingly  dominate  the  response. 
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Figure  7.  Force  versus  displacement,  effects  of  Young’s  modulus  E  for  R  =  25  jum :  (a)  d  =  5  /an,  (b)  d  =  0.05  //m . 


Figure  8.  Indentation  force  P  versus  displacement  S,  effects  of  Poisson’s  ratio  V  for  R  =  25  //m:  (a)  d  =  5  //m, 
(b)  d  =  0.05  jum. 


Figure  9  shows  the  effect  of  yield  strength  on  the  P-8  response.  Two  data  sets  are  presented, 
corresponding  to  the  actual  titanium  alloy  of  figure  3  (black  curve),  and  a  hypothetical  material 
with  the  same  elastic  properties  as  the  former  but  with  a  flow  stress  reduced  by  50%  (red  curve). 
We  compare  figures  7(a)  and  9,  a  50%  reduction  in  flow  stress  seems  to  exert  much  more 
significant  effect  on  the  indentation  force  measured  during  the  loading  phase  of  the  simulation 
than  does  a  50%  reduction  in  Young’s  modulus.  On  the  other  hand,  the  substrate  material’s  flow 
stress  apparently  imposes  a  negligible  effect  upon  the  force  measured  during  the  unloading 
portion  of  the  cycle. 

Notice  also  from  figures  7  through  9  the  fairly  jagged  nature  of  the  load  versus  displacement 
curves,  particularly  at  shallow  indentation  depths.  This  phenomenon  is  not  attributable  to  the 
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discreteness  of  the  data,  since  each  curve  in  each  figure  is  comprised  of  hundreds  of  data  points, 
Instead,  it  is  thought  to  be  a  result  of  the  discreteness  of  the  indenter  geometry,  constructed  of  a 
series  of  line  segments  in  the  finite  element  model.  The  jaggedness  was  not  alleviated  upon 
mesh  refinement,  the  use  of  “full  integration”  as  opposed  to  “reduced  integration”  elements,  or 
variation  of  the  hourglass  stiffness  parameter  used  in  reduced  integration  elements  to  suppress 
excessive  distortion.  A  second  cause  could  be  the  presence  of  numerical  instabilities  in  the 
implicit  calculations  (e.g.,  loss  of  positive-definite  tangent  stiffness  because  of  relatively  weak 
strain  hardening  prescribed  in  the  constitutive  model). 


£>[//m] 

Figure  9.  Indentation  force  P  versus  displacement  &  effects  of  flow 
stress  for  R  =  25  jLttn  and  d  =  5  fum  and. 


4.  Experiments  and  Simulations:  A  Comparison  of  Results 


4.1  Ti-6A1-4V  Alloy 

Simulated  and  experimentally  obtained  results  are  compared  for  the  Ti-6A1-4V  alloy  over  a  range 
of  indentation  depths  in  figure  10.  Physical  experiments  were  conducted  via  the  instrumented 
indentation  technique  of  VanLandingham  (2003),  featuring  a  diamond-tipped  spherical  indenter 
(radius  R  =  25  //m)  undergoing  small  dynamic  oscillations  superposed  upon  the  macroscopic  load- 
unload  cycle.  From  figure  10,  the  P-fi curves  agree  well  during  the  loading  phase;  however,  the 
simulations  exhibit  a  noticeably  stiffer  response  upon  unloading  than  do  the  experiments.  This 
discrepancy  may  be  attributed  to  elastic  and/or  plastic  anisotropy:  the  constitutive  model  used  in 
the  simulations  is  isotropic,  while  Ti-6A1-4V  is  a  hexagonal  microstructure  exhibiting  an  aniso¬ 
tropic  tangent  stiffness  at  the  single  crystal  level  (cf.  Schoenfeld  &  Kad,  2002).  Furthennore,  local 
effects  of  microstructural  heterogeneity  at  the  scale  of  the  contact  region  (e.g.,  grain  and  phase 
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boundaries,  dislocation  “pile-ups”  or  substructure)  are  not  explicitly  captured  by  the  classical 
plasticity  model  nor  are  directional  strain  history  effects  (e.g.,  kinematic  hardening  because  of 
stress  directionality  under  non-monotonic,  multi-axial  deformation  modes).  The  former  effects 
(i.e.,  heterogeneity  of  microstructure)  are  thought  to  be  of  relatively  low  importance  here,  even 
though  the  grain  size  is  on  the  order  of  the  indenter  size,  since  experimental  repeatability  was 
deemed  satisfactory.  Finally,  the  experimental  input  data  for  Ti-6A1-4V  were  acquired  from  a 
uniaxial  tension  test  (figure  3),  while  the  indentation  experiment  produces  a  primarily  compressive 
stress  state.  Thus,  discrepancies  could  be  caused  by  tension-compression  asymmetry  of  the 
material  response.  Note  that  the  overshoot  of  the  experimental  curves  by  the  simulation  curves 
during  the  loading  phase  does  not  indicate  a  deficiency  in  predictive  capability  of  the  model  but 
merely  a  differently  specified  final  indentation  depth  d  than  that  attained  in  the  experiment. 


0.0  0.5  1.0  1.5  2.0 

8\ju  m] 


Simulation  parameters 

E  =  112  GPa 
v=  0.31 

cr  (0)  =  1000  MPa 

y 


Figure  10.  Indentation  force  P  versus  displacement  &  Ti-6A1-4V,  experiment  and 
simulation. 


Figure  1 1  (red  curve)  shows  the  stress-strain  prediction  reproduced  from  results  of  a  simulation, 
conducted  with  the  material  properties  given  in  figure  3,  with  these  input  parameters  corre¬ 
sponding  to  the  black  curve  in  figure  1 1 .  The  prediction  (red  curve)  is  based  on  the  analytical 
model  of  Herbert  et  al.  (2001),  wherein  the  equivalent  true  stress  a  and  strain  £  are  defined  as 

,  P  .  P  ,  0.2 a.  0-2 fidjEdl  m 

^  jt(2  djt-d/y  R  R  w 

in  which  a  is  the  radius  of  contact  and 

d<=s~%  (10) 

is  the  estimated  contact  depth  for  a  spherical  indenter.  The  elastic  contact  stiffness  S  is  given  by 

S  =  2aEr,  (11) 
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and  the  reduced  elastic  modulus  is 


i-y, , 

E.  E. 


(12) 


in  which  E,-  and  Es  are  Young’s  moduli  of  the  indenter  and  substrate,  respectively,  and  v,  and  Vs 
are  Poisson’s  ratios  of  the  indenter  and  substrate,  respectively.  Here,  for  a  rigid  indenter,  Er  — > 
Es/(  1  -  vs).  Given  P  as  a  function  of  8,  equations  9  through  1 1  are  solved  iteratively  for  <7  and 
£  .  Agreement  between  original  (black)  and  reproduced  (red)  curves  is  considered  here  to  be 
marginal  at  best,  with  the  initial  yield  stress  drastically  under-predicted  at  small  strains  (i.e.,  for 
£  <  0.05  ).  Herbert  et  al.  (2001)  reported  a  similar  difficulty  in  matching  uniaxial  and  indenta¬ 
tion  data  for  6061-T6  aluminum.  Also  shown  in  figure  1 1  is  Hertz’s  elastic  solution,  in  which 
the  force  is  calculated  by 


P  = 


4  Era3 
3  R 


(13) 


and  is  then  substituted  into  equations  9  through  1 1 .  The  elastic  response  of  the  material  is 
reproduced  quite  well  by  the  theory;  of  course,  equation  13  is  invalid  upon  yield,  which  occurs 
here  at  an  applied  strain  of  £  —  0.0089  . 


Figure  11.  Experimental,  simulated,  and  analytical  elastic 
uniaxial  stress-strain  responses. 


4.2  W-Sn  Alloys 

The  indentation  responses  of  two  different  W-Sn  metallic  materials  were  examined,  labeled  here 
as  “70K32”  and  “K4754”.  As  shown  in  figure  12,  Young’s  modulus  E  and/or  yield  strength  oy 
were  varied  parametrically  in  the  simulations.  Elastic-perfectly  plastic  behavior  was  assumed  in 
all  cases,  as  was  a  fixed  Poisson’s  ratio  of  v=  0.4.  Insufficient  experimental  data  were  available 
to  justify  implementation  of  a  more  elaborate  constitutive  model,  accounting  for  strain  hardening 
or  viscous  effects.  Ranges  of  yield  strengths  and  elastic  moduli  were  chosen  on  the  basis  of 


11 


values  determined  by  uniaxial  tension  and  instrumented  dynamic  indentation  experiments.  For 
the  alloy  70K32  of  figure  12(a),  it  appears  that  the  assumption  of  constant  yield  stress  (i.e., 
perfect  plasticity)  may  be  inadequate,  as  the  slope  of  the  experimental  curve  during  the  loading 
phase  of  the  indentation  cycle  increases  with  depth  <7,  while  the  simulated  results  predict  a  P-o 
response  with  little  curvature.  On  the  other  hand,  the  load  displacement  response  for  alloy 
K4754  of  figure  12(b)  appears  to  be  well  characterized  by  an  elastic-perfectly  plastic  model.  For 
both  W-Sn  materials,  the  elastic  modulus  E  appears  to  exert  little  influence  on  the  response 
relative  to  the  yield  strength  <7V  at  depths  greater  than  several  hundred  nanometers. 


-  Experiment 

-  E  =  50  GPa,  cry  =  43  MPa 

E  =  35  GPa,  cry  =  43  MPa 
E  =  20  GPa,  (jy  =  43  MPa 

-  E  =  35  GPa,  cry  =  47  MPa 

E  =  35  GPa,  (Ty  =  39  MPa 


£[nm] 

(a) 


Experiment 

E  =  50  GPa,  ay  =  40.8  MPa 
E  =  35  GPa,  rry  =  40.8  MPa 
E  =  20  GPa,  r7y  =  40.8  MPa 


£>[nm] 

(b) 


Figure  12.  Indentation  force  P  versus  displacement  8 for  (a)  W-Sn  alloy  70K32  and  (b)  W-Sn  alloy  K4754. 


4.3  W-Nylon  Blends 

The  indentation  responses  of  three  distinct  W-Nylon  composites  were  studied,  labeled  as 
“RVY4B,”  “026QX,”  and  “1 AA63”.  These  are  high  density  metal-polymer  mixtures. 
Homogeneous  elastic-perfectly  plastic  behavior  was  assumed  in  the  simulations  of  these 
materials,  again  with  a  fixed  Poisson’s  ratio  of  V=  0.4.  No  composite  theory  (i.e.,  properties 
weighted  by  volume  fractions  of  W  and  nylon)  was  used.  As  shown  in  figure  13(a),  Young’s 
modulus  E  and/or  yield  strength  o\,  were  varied  parametrically  in  the  simulations  of  material 
RW4B.  Ranges  of  yield  strengths  and  elastic  moduli  were  again  chosen  based  upon  values 
predicted  by  uniaxial  tension  tests  and  instrumented  dynamic  indentation  experiments.  Clearly, 
from  figure  13(a),  experiment  and  simulation  differ  dramatically  for  all  values  of  simulation 
parameters.  It  is  hypothesized  that  the  elastic-perfectly  plastic  model  used  in  these  simulations 
is  inadequate  because  of  its  inability  to  capture  strain  rate  effects  (i.e.,  viscoelasticity,  time 


12 


dependence,  and  hysteresis)  engendered  by  nylon  components  of  the  microstructure.  Single 
simulations  were  conducted  for  W-Nylon  materials  026QX  and  1AA63.  Because  of  the 
inability  of  the  constitutive  model  to  accurately  represent  the  behavior  of  these  systems,  ranges 
of  E  and  <TV  were  not  explored  for  these  materials. 


c>[nm] 

(a) 


J[nm]  £[nm] 


(b) 


(c) 


Figure  13.  Indentation  force  P  versus  displacement  5 for  W-nylon  materials  (a)  RVV4B,  (b)  026QX,  and  (c)  1AA63. 


13 


5.  Conclusions  and  Recommendations 


An  axisymmetric  model  of  spherical  indentation  has  been  developed  and  invoked  within  an 
implicit  Lagrangian  finite  element  setting.  The  model  was  used  to  study  the  load  displacement 
response  of  elastic-plastic  substrates  subjected  to  quasi-static  contact  with  a  rigid  indenter.  The 
implemented  constitutive  model  assumed  rate-independent,  isotropic  elastic-plastic  behavior  of 
the  substrate  materials.  Materials  studied  included  Ti-6A1-4V  alloy,  W-Sn  alloys,  and  W-Nylon 
composites.  Numerical  predictions  were  compared  with  results  obtained  from  instrumented 
indentation  experiments.  Key  findings  regarding  specific  materials  include 

•  Indentation  experiment  and  simulation  are  capable  of  close  agreement  in  tenns  of  load 
displacement  response  when  the  constitutive  model  used  in  the  simulation  is  accurate  with 
regard  to  phenomenology  and  material  parameters,  as  in  the  case  here  of  a  Ti-6A1-4V  alloy. 
For  this  material,  the  loading  portion  of  the  experimental  response  was  duplicated 
accurately  via  finite  element  simulations,  while  the  elastic  unloading  portion  was  not  as 
accurately  reproduced.  Assumptions  of  frictionless  contact  and  a  rigid  indenter  did  not 
inhibit  simulation  accuracy. 

•  Classical  rate-independent,  isotropic  elastic-plastic  constitutive  models  appear  adequate  for 
modeling  the  response  of  the  featured  W-Sn  materials.  The  alloy  70K32  appears  to 
undergo  strain  hardening  during  indentation,  whereas  the  alloy  K4754  appears  to  exhibit 
perfectly  plastic  behavior. 

•  Classical  rate-independent,  isotropic  elastic-plastic  constitutive  models  appear  inadequate 
for  modeling  the  response  of  the  featured  W-Nylon  materials.  These  materials  seem  to 
exhibit  viscoelastic  behavior,  necessitating  an  accurate  treatment  of  time-dependent  finite 
elastic  defonnations. 

Table  1  lists  material  parameters  provided  by  experiment  (Young’s  modulus  E  from  indentation 
tests;  initial  yield  stress  <Jy  from  uniaxial  tension  tests)  and  choices  of  these  parameters  that 
enabled  simulations  to  best  match  the  experimental  indentation  data.  Properties  prescribed  with 
relative  confidence  are  limited  to  yield  strengths  for  the  W-Sn  materials,  highlighted  in  yellow. 
The  present  set  of  simulations  did  not  permit  accurate  prescription  of  properties  for  the  W-Nylon 
materials;  however,  the  simulations  did  permit  us  to  substantiate  the  inadequacy  of  rate- 
independent  metal  plasticity  constitutive  models  for  describing  their  mechanical  behavior.  Also 
shown  for  comparison  purposes  are  the  ideal  properties  used  for  the  Ti-6A1-4V  material  used  to 
validate  the  computational  approach. 
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Table  1.  Material  parameters  obtained  experimentally  and  verified  by  indentation  simulation. 


Material 

Identifier 

Composition 

E  (GPa) 

<jv  (MPa) 

Experiment 

Simulation 

Experiment 

Simulation 

70K32 

Tungsten-Tin 

35  ±  15 

35 

43  ±3.6 

46.6 

K4754 

Tungsten-Tin 

35  ±  15 

35 

40.8  ±7.4 

40.8 

RVV4B 

Tungsten-Nylon 

4.2  ±6 

— 

28.5  ±3.7 

— 

026QX 

Tungsten-Nylon 

2.3  ±3 

— 

29.2  ±4.4 

— 

1AA63 

Tungsten-Nylon 

1.6  ±2 

— 

36.7  ±4.3 

— 

Ti-6A1-4V 

Titanium- Aluminum-  Vanadium 

112 

112 

1000 

1000 

In  summary,  for  homogeneous  isotropic  elastic-plastic  substrate  materials  exhibiting  small 
elastic  strains  (e.g.,  most  conventional  polycrystalline  metals),  coupled  experimental  and 
numerical  indentation  methods 

•  Appear  to  enable  verification  of  Young’s  modulus  measurements  for  elastic  (shallow) 
indentation. 

•  Appear  to  enable  verification  of  yield  strength  measurements  for  plastic  (deep)  indentation. 

•  Do  not  appear  to  enable  verification  of  Poisson’s  ratio,  unless  Young’s  modulus  is  known 
precisely. 

•  Appear  to  enable  one  to  detect  the  presence  of  and  estimate  the  degree  of  inelastic  strain 
hardening. 

In  future  studies  of  metals,  it  may  be  worthwhile  to  implement  crystal-plasticity-based  formu¬ 
lations  (cf.  Schoenfeld  &  Kad,  2002)  accounting  for  anisotropic  elastic  and  inelastic  responses 
engendered  by  local  crystallographic  orientation,  since  the  contact  region  for  micron-sized 
diameters  will  typically  be  confined  to  surface  areas  spanning  one  or  a  few  crystals  in  typical 
engineering  metals  with  grain  sizes  on  the  order  of  tens  of  microns.  Analogously,  more  descrip¬ 
tive  constitutive  theories  for  polymers  may  also  offer  closer  agreement  between  simulation  and 
experiment  (e.g.,  models  capturing  finite  anisotropic  viscoelastic  effects)  (cf.  Kaliske,  2000). 

Note  that  at  even  finer  scales  of  spatial  resolution,  microstructural  heterogeneities  in  the  form  of 
dislocations,  grain  and  phase  boundaries,  and  point  defects  may  exert  considerable  influence  on 
the  response,  especially  for  nanoscopic  indenter  dimensions.  In  such  situations,  continuum 
mechanics  treatments  may  be  merged  with  atomic  physics  in  order  to  resolve  the  effects  of 
individual  lattice  defects,  for  example,  in  coupled  atomistic-continuum  schemes  such  as  the  quasi¬ 
continuum  method  (Tadmor  et  ah,  1996;  Shenoy  et  ah,  1999)  or  atomistically  infonned  fully 
finite  element  methods  (Zhu  et  ah,  2004). 
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Appendix  A.  ABAQUS  Input  Deck 


An  ABAQUS  input  deck  is  given  below  for  quasi-static  indentation  and  release  of  a  slab  of  Ti- 
6A1-4V  material,  with  an  indenter  radius  of  25  //m ,  to  an  indentation  depth  of  5  f.im  .  The 

software  used  for  the  present  work  was  version  6.3-5  (2002). 


*HEADING 

STATIC  INDENTATION  WITH  A  HEMISPHERICAL  PUNCH 

■k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

**  Measure  the  punch  reaction  force  during  the  following  2  steps: 

**  Step  1:  Displace  punch  1R  units  downwards. 

**  Step  2:  Return  punch  to  original  position. 

**  (Force  Units:  N,  MPa,  mm,  sec) 

**  R  is  the  unit  radius  of  the  punch. 

**  Mesh  is  10R  x  10R,  1/5  suggested  by  Hill  et  al.,Pr.  Roy.  Soc.,1989 
**  4900  elements,  5041  nodes  for  target  material 

■k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k'k 

*PREPRINT, MODEL=YES, HISTORY=NO 

*RESTART, WRITE 

*NODE 

1,  0 . ,  250  . 

71,  250. ,250. 

4971,  0.,0. 

5041,250 . ,  0  . 

*NSET, NSET=N1 

U 

*NSET, NSET=N7 1 
71, 

*NSET, NSET=N4  97 1 
4971, 

*NSET, NSET=N50  4 1 
5041, 

*NFILL, NSET=TOP, BIAS=0 . 95 
N1 , N7 1 , 70, 1 

*NFILL,  NSET=BOT, BIAS  =  0 . 95 
N4971, N5041, 70, 1 
*NFILL, NSET=ALLN, BIAS=0 . 95 
TOP, BOT, 70, 71 

*NSET, NSET=CENTER, GENERATE 

1.4971.71 

*NSET, NSET=OUTSIDE, GENERATE 

71.5041.71 
*ELEMENT, TYPE=CAX4R 
1, 72, 73,  2,  1 

*ELGEN, ELSET=SUBSTRATE 

1.70. 1. 1.70.71.70 
*ELSET, ELSET=CENT, GENERATE 

1. 4901. 70 

*ELSET, ELSET=ETOP, GENERATE 
1,70,  1 
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*NODE, NSET=PUNCH 
100000,0. ,275. 

*SOLID  SECTION, ELSET=SUBSTRATE, MATERIAL=TIAL 
* HOURGLASS  STIFFNESS 
213 . 4 


* MATERIAL, NAME=TIAL 

*ELASTIC 

112000,  0.31 

*PLASTIC,  HARDEN ING= ISOTROPIC 


1000, 

o 

o 

1075. 

9739, 

0 

.  0201 

1099. 

8  646, 

0 

.  0324 

1125 . 

0435, 

0 

.  0526 

1150  . 

9742, 

0 

.  0793 

1175  . 

6331, 

0 

.  1024 

1200  . 

0284, 

0 

.  1344 

1225  . 

2135, 

0 

.  1738 

1250  . 

0531, 

0 

.2272 

1275  . 

4083, 

0 

.2908 

1290  . 

6458, 

0 

.3369 

* BOUNDARY 
BOT, 1,2 
CENTER,  1 
100000,  1 
100000, 6 
OUTSIDE, 1 


*SURFACE,  TYPE=SEGMENTS, NAME=IMP ACTOR 
START, 25 . , 275 . 

CIRCL, 0 . , 250 . , 0 . , 275 . 

* SURFACE,  TYPE=ELEMENT, NAME=TARGET 
ETOP, S3 

*RIGID  BODY, ANALYTICALSURFACE=IMP ACTOR, REFNODE=100000 
*CONTACT  PAIR, INTERACT I ON=ROUGH 
TARGET, IMP ACTOR 

*SURFACE  INTERACTION, NAME=ROUGH 
*FRICTION 

0 .0, 


*  STEP ,  NLGEOM,  INC=2  0  0 , AMPLI TUDE=RAMP , UNS YMM=YES 
Displace  punch  XR  units  downwards. 

*STATIC 
. 0015, 1 .  ,  ,  .  05 
* BOUNDARY 
100000,2,2,-5 
*PRINT, CONTACT=YES 

*CONTACT  CONTROLS, FRICTION  ONSET=DELAY 
*CONTACT  PRINT, SLAVE=TARGET 

*  ^CONTACT  FILE, SLAVE=TARGET, FREQUENCY=10 
*EL  PRINT, FREQUENCY=50 , ELSET=CENT 

S, 

E, 

*NODE  PRINT, FREQUENCY=1 
U, 

RF, 
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*OUTPUT, FIELD, FREQUENCY=1 
*NODE  OUTPUT 


U,  V,  A,  RF 
*ELEMENT  OUTPUT 
S, LE, NE 

*CONTACT  OUTPUT, VARIABLE=PRESELECT,  NSET=TOP 
*OUTPUT, HISTORY, FREQUENCY=1 
*NODE  OUTPUT, NSET=PUNCH 
U,  V,  A,  RF 


*ELEMENT  OUTPUT, ELSET=CENT 
S, LE, NE 
*END  STEP 

■k-k-k-k-k-k-k'k-k-k-k-k-k-k'k-k-k-k-k-k-k-k-k-k-k-k-k-k-k-k-k'k-k-k-k-k'k-k-k'k-k-k-k-k-k-k'k-k-k-k-k-k 


*  STEP ,  NLGEOM,  INC=2  0  0 , AMPLI TUDE=RAMP , UNS YMM=YES 
Return  punch  to  original  position. 

*STATIC 

.0015, 1 .  ,  ,  .015 
* BOUNDARY, OP=MOD 
100000,2,2,0.0 
*OUTPUT, FIELD, FREQUENCY=1 
*NODE  OUTPUT 
U,  V,  A,  RF 
*ELEMENT  OUTPUT 
S, LE, NE 

*CONTACT  OUTPUT, VARIABLE=PRESELECT, NSET=TOP 
*OUTPUT, HISTORY, FREQUENCY=1 
*NODE  OUTPUT, NSET=PUNCH 
U,  V,  A,  RF 

*ELEMENT  OUTPUT, ELSET=CENT 
S, LE, NE 

*ENERGY  OUTPUT, VARIABLE=PRESELECT 
*END  STEP 
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